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ABSTRACT 

We present a parameter study of the possibility of tidally triggered disk instability. 
Using a restricted iV-body model which allows for a survey of an extended parameter 
space, we show that a passing dwarf star with a mass between 0.1 and 1 Mq can 
probably induce gravitational instabilities in the pre-planetary solar disk for prograde 
passages with minimum separations below 80-170 AU for isothermal or adiabatic 
disks. Inclined and retrograde encounters lead to similar results but require slightly 
closer passag es. Such encounter distances are quite likel y in young moderately massive 
star clusters l|Scallv fc Clarkel200lUBonnell et all2001|) . The induced gravitational in- 
stabilities may lead to enhanced planetesimal formation in the outer regions of the 
protoplanetary disk, and could therefore be relevant for the existence of Uran us and 
Neptune, whose formation timescale of about 100 Myr ijWuchterl et alJ l2"000) is in- 
co nsistent with the dis k lifetimes of about a few Myr according to observational data 
bv lHaisch et alJ l)2001|) . The relatively small gas/solid ratio in Uranus and Neptune 
can be matched if the perturbing fly-by occurred after early gas depletion of the solar 
system, i.e. when the solar system wa s older than a bout 5 Myr. 

We also confirm earlier results bv lHellerl l|1993t) that the observed 7 degree tilt of 
the solar equatorial plane relative to the ecliptic plane could be the consequence of 
such a close encounter. 

Key words: stellar dynamics - methods: N-body simulations - Kuiper Belt - minor 
planets, asteroids - Solar system: formation - open clusters and associations: general 



1 INTRODUCTION 



1.1 Limits of the coagulation model 

In the classical model of planet formation dust grains col- 
lide and stick together forming larger grains and clumps. 
These clumps also collide with each other forming progres- 
sively larger planetesimal bodies. At some point the mass of 
a planetesimal provides sufficient gravitation to attract sur- 
rounding particles — the accretion process begins. While the 
coagulatio n process takes several Myr fo r Jupiter and Sat- 
urn llThommes. Duncan_fc_Le vison 200S ), it t akes at least 
100 Mvr for Neptune faollenbach et alJ l200(f> and proba- 
bly even longer for objects in the Edgeworth-Kuiper Belt. 
The timescale T coag of the coagulation phase is mainly de- 
termined by the density g of the surrounding protoplanetary 
dust. It further depends on the epicyclic frequency k since 
the epicyclic motion is an important factor for the collision 
rate of planetesimals. k, is given by 



r |-4fi 

dr 



(1) 



where fi is the orbital angular frequency and r the ra- 
dius of the particle orbit, i.e. the distance from the central 
star. The dust density g is approximately proportional to 
the surface density S , which is expected to fit a r~ 3 ^ 2 law 
JWeidenschillinJl977llMaver et all2004h . Thus for low mass 
disks 
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(2) 



where f^Kepier is the Keplerian frequency. The relation- 
ship between T coag and r in a gaseous disk containing 
fine dust with negligible gravitatio nal interaction between 
the dust grains can be written as <ISafronovlll96l see also 
iRice fc Armitagdl200 j) 



(3) 
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The relation may be somewhat less steep due to the smaller 
end masses of the outer planets and gravitational focus- 
ing for large r dust part icles. This is i n agre ement with 
iPollack et al] (Il996l) and iKokubo fc Idal l200d) who found 
Tform oc r v with eta — 2 and r\ « 2.6, respectively. Assum- 
ing a typical formation timescale for Jupiter's core of 3 Myr 
(''jup = 5.2 AU) this relation leads to prohibitively long for- 
mation timescales up to 600 Myr. Even the conservative r 2 
law yields timescales of 40 Myr for Uranus and 100 Myr for 
Neptune. Thus we find good agreement of the formation- 
time scaling given by Eq. [3]witli the more detailed investi- 
gations discussed at the beginning of this section. 

Recent work bv IRafikovl i2004h predicts shorter forma- 
tion times in the order of ~ 10 7 yr for the core of an ice giant 
at 30 AU from the Sun. Although this is only one-sixth of the 
conventional value, it is also near the disk lifetime discussed 
below and even beyond, if the accretion time of the envelope 
and t he precursory depletion of large plan etesimals (>0.1- 
1 km ; lRafikovll20 04) is taken into account. iGoldreich et all 
(2004) also predict similarly short formation timescales but 
require a six times higher surface density ( compared to the 
minim um mass solar n ebula, MMSN; see IWeidensc hilling 
119771 and lHavasrjll98l in the region of Uranus and Nep- 
tune to allow them to collect most of their final mass until 
the end of oligarchy. This would also imply a significantly 
higher total disk mass. 

The lifetime of the protoplanetary disk is most proba- 
bly less than 10 Myr . Base d on JHKL-excess observations 
lHaisch. Lada fc Ladal ll200lf) found lifetimes less than 6 Myr 
before the disk material is dispersed by the stellar wind, pho- 
toevaporation and gravitational scattering by protoplanets. 
The masses of observed disks range from .003 to 0.3 Mq, 
with the bulk below 0.03 Mq llNattall2004h . This is compa- 
rable to the MMSN which is between 0.01 and 0.07 Mq. 
Thus, unless an unusual high disk mass or surface density 
profile is assumed, there is a serious problem to explain the 
existence of Uranus and Neptune. There have been several 
theoretical attempts to solve this contradiction: 



(i) The coagulation and accretion process could be much 
faster than predicted by current models. But it seems un- 
likely that the models overestimate the formation times by 
more than an order of ma gnitude. Even the fast accretion 
scenario by IRafikovl d2004l) can hardly solve this inconsis- 
tency. 

(ii) Protoplanetary disks could survive much longer than 
expec ted from observations like those by lHaisch et alJ 
(2001), since only the dust component of a disk can be de- 
tected by current observational methods, but no solid plan- 
etesimals. But as Uranus and Neptun e mainly consist of soli d 
material and only about 10% of gas dWuchterl et alj feoOO'). 
this does not seem to be a satisfactory attempt either. 

(iii) Uranus and Neptune may have formed much closer 
to the Sun in the birth region of Jupiter and Saturn and 
ejected from their initi al orbits after reaching a certain mass 
dThommes et al.ll2002T) . With this scenario, however, it is 
difficult to understand the almost circular actual orbits of 
Neptune and Uranus. 



1.2 Gravitational instabilities 

Fragmentation due to gravitational instabilities (GIs) is an- 
other possible solution to this problem. The dynamical free- 
fall timescale of a GI in a cloud with a density excess g exc 
(relative to the ambient medium) can be estimated as 



Iff = 



3tt 



32G£> C 



(4) 



dBinnev fc Tremaind fl987f) . where G is the gravitational 
constant. For an over-density g cxc = 10~ 12 g cm -3 (which 
would result by doubling the ambient density locally in 
the disk at 30 AU, the today's orbital radius of Neptune) 
Tff w 100 yr, which is six orders of magnitude shorter than 
the classical formation timescale and can therefore be ne- 
glected compared to the subsequent accretion phase. 

There are two major ways for the development of GIs. A 
very massive disk of more than 0.1 Mq becomes unstable to 
axisymmetric perturbati ons according t o the disk stability 
criterion Q developed bv lToomrel dl964l) . 
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where c s is the speed of sound in the gas, k the epicyclic 
frequency and E the surface density. 

A disk can become unstable if it continues to ac- 
crete material from the protostellar cloud, so that Q be- 
comes less than 1. For non-axisymmetric perturbations in 
aseous disks, the limit is higher at Q ~ \/3 « 1.7 
Polvachenko. Polvachenko fc Strel'Nikovlll997t) . 

Another way towards local instabilities is via tidal per- 
turbations. With this contribution we argue that also ini- 
tially stable disks can become locally unstable if an external 
perturbation on the surface density of the disk causes Q to 
fall below 1 locally, so that fragmentation only occurs in 
such regions. For example and as shown below, a star with 
mass m = 0.5 Mq on a hyperbolic orbit with eccentricity 
e = 1.5 and pericentre distance tq = 120 AU may cause a 
local Toomre instability within 500 years. Such encounters 
are expected to occur with a realistic proba bility within the 
first 6 Myr in the Orion Nebula Cluster l|Scallv fc Clarke! 
120011: iBonnell et al.ll200lT) . 

The tilt of the solar rotational axis aga inst the nor- 
mal direction of t he ecliptic plane of about 7° l|Eggers et alJ 
ll997tlHellerlll993Fl is int eresting in thi s respect. We confirm 
the results presented bv iHellen l)l993F) . who used smoothed 
particle hydrodynamics (SPH) instead of our restricted TV- 
body calculation, that such an inclined encounter leads to 
the observed tilt. 

In this paper we study if gravitational instabilities, as 
a possible originator of planet formation, can be triggered 
due to stellar fly-bys in a parameter survey of the stellar 
mass, encounter distance and eccentricity, m, vq and e, re- 
spectively. The aim of this survey is to identify those regions 
in the m-ro-e parameter space where growing local instabil- 
ities can be triggered. Due to the extended parameter space 
such a study can only be performed by fast numerical meth- 
ods. Though detailed hydrodynamical calculations would be 
preferable for such models, they are too CPU expensive for 
scanning a large set of parameter space. Therefore, we ap- 
plied a simple numerical ansatz based on the fast restri cted 
N-body method as was already used bv llda et al.ll|2000|) . In 
the future we will use our present survey to perform hydro- 
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dynamical calculations for a few interesting m, ro and e in 
order to investigate the scenario of tidally triggered GIs in 
more detail. 

In §|5|we describe our numerical model. §[3]describes the 

general effects of star-star encounters on circumstellar disks 
for an individual example. §3] gives a plan of the parameter 
study and shows the results for coplanar and non-coplanar 
encounters. Orbits of protoplanet candidates and the tilt of 
the circumstellar disk due to inclined fly-bys are also de- 
scribed. In §|K|the results are discussed and the conclusions 
are presented in § [(J 



1.3 Encounter probability 

Tidally-induced gravitational instabilities can only be a vi- 
able mechanism if encounters are of sufficient likelihood. To 
estimate the encounter probability we need to consider a 
model of a young cluster. 

Assuming a Plummer star-cluster model with given 
half-mass radius P0.5 = -Rpi/\/2 2 / 3 — 1 ~ 1.305 Rp\ (Rp\ is 
the Plummer-radius), mass M c \ and number of stars N c i, 
we can calculate the crossing time T cr from the character- 
istic velocity dispersion and compute the expected number 
of encounters n cnc within a given impact parameter b and 
therefore the mean time between encounters. This gives us a 
useful guide-line as to how likely the interesting encounters 
may be. 

The relationship between 6, the fly-by eccentricity e and 
the pericentre distance ro is 



ro = 



1 



e + 1 



(6) 



The one-dimensional characteristic velocity dispersion ctid 
in the Plummer model, derived from the virial theorem, is 
given by 
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and the crossing time is 
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This leads to the mean time between encounters 
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n cnc N c \ 

The probability P onc for at least one event within the disk 
lifetime Td can be derived from the Poisson-probability for 
x encounters, 



P ell c(x) = — exp (— v) , where v 
x>. 

For at least one event, 

Pcnc(x >!) = !- Pcnc(O) = 1 - exp 



-P^ (10) 
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For a young cluster with M c \ = 500 Mq, R0.5 ~ 0.3 pc 
and with a mean stellar mass of 0.5 Mq we 
get 10 % < Pone < 32 % within T d = 6 Myr for 
60 AU < r < 200 AU, as shown in Fig □ This means 
that such close encounters are indeed likely in young stellar 
clusters. The chances are even greater, up to 50%, for a 
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Figure 1. Probability P en c of close encounters within = 6 Myr 
as a function of the minimum encounter distance ro (Eq. |S| for 
three different Plummer clusters: A young open cluster of high 
density (solid curve), a more expanded cluster of the same mass 
(thick dotted curve) and a young ONC-type cluster (thin dotted 
curve, iKroupa et al 120011) . Pcnc increases with increasing cluster 
mass M c \ and decreasing half-mass radius. 



cluster of 5000 Mq and Pp. 5 = 0.5 pc, similar to the Orio n 
Nebula Cluster (ONC ) iKroupa. Aarseth fc HurTevll200ll) . 
IScallv fc Clarke! (l200ll) propose an encounter probability of 
only 4 % for encounters closer than 100 AU in the ONC at 
its present density and age. This corresponds to a proba- 
bility of 9 % for ro < 150 AU, and is still a reasonably high 
likelihood. We note that the probabilities for encounters in 
a younger and denser ONC would be accordingly higher. 



2 THE MODEL 

We investigate a large region of parameter space in order 
to study in which fraction of it tidally induced GI may be 
relevant. Our parameter space is defined by the distance 
of closest approach ro, the eccentricity of the relative or- 
bit e and the mass of the perturbing star m. While an ex- 
act description of gravitational instabilities requires a 3d- 
stellar-hydrodynamical calculation, we apply a simpler, but 
numerically much faster restricted N-body approach. This 
ansatz is valid as long as the self-gravity of the disk, the 
pressure forces or the inelastic collisions between the con- 
stituents of the disk can be neglected. For example, simu- 
lations of encounters betwee n galaxies have demonst rated 
the validity of this approach (iToomre fc Toomreill972}l . Be- 
cause 3d-hydrodynamical calculations are by far too CPU 
intensive at the moment, such simplified dynamical mod- 
elling is the only way to scan the parameter space in detail. 
Of course, we cannot follow the dynamical evolution of sub- 
structures when their self-gravity becomes important or even 
dominant. However, we can identify regions which are prone 
to such instabilities. 

In our computations the disk is represented by test par- 
ticles on initially circular and coplanar orbits. There is no ex- 
plicit interaction (self- gravity, coagulation, pressure forces) 
between the particles themselves and also no gravitational 
feedback of the disk to the Sun or the perturbing star. We 
assume both stars to move on hyperbolic orbits. This should 
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be typical for young star clusters. Our numerical method re- 
quires low computational effort. Therefore, it allows for the 
requested in-depth parameter study of the main orbital pa- 
rameters m, e and ro providing an estimate of the order 
of magnitude of the perturbation and its probable effects on 
disk stability. We expect this approach to give useful insights 
into those regions of the encounter-parameter space where 
the physics will be of interest for triggered GI, because the 
primary mechanism driving the fluctuations in disk density 
within the first 500 years comes from the changes of the 
overall potential due to the relative orbit of the two stars. 



2.1 Disk properties 

In our model we suppose an initially flat disk with an outer 
radius r max = 100 AU, an inner border r m i n = 3AU and a 
total mass M<j = 0.073 Mq, which is close to the Minimum 
Mass Solar Nebula. The surface density E is then given by 



E(r)= 6300 gem" 2 - 



-3/2 



(12) 



For later analysis (§ 14. H the disk is as sumed to co ntain a 
dust fraction of the order of 1 per cent jNattall2004fo . Thus 
the surface density of the dust component is approximately 
l/100th of the total E, E du st(r) « 0.01E(r). The density 
profile is needed for the calculation of the local Toomre pa- 
rameter Q which we use as an indicator of local gravitational 
instability. As shown in Eq. Q>] Q depends on the speed of 
sound in the gas fraction of the disk, which mainly depends 
on the radial temperat ure profile. Following a recommenda- 
tion bv lRafikovl ||2004) we use 



1200 ms" 



(w) 



Thus, the initial Q is approximately given by 
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The changes of Q also depend on the adopted polytropic 
index 7, which is given with respect to the cooling efficiency 
of the disk. Large cooling rates prevent tidally compressed 
clouds from significant heating, i.e. we have an isothermal 
gas, 7 = 1. If the cooling is negligible the gas is compressed 
adiabatically. For a di-atomic gas 7 = 1.4 which is a good 
approximation for the protoplanetary gas since it mainly 
consists of molecular hydrogen. Since p tx p 7 and c s tx \fp~j~p 
and E tx p, the dependence of Q on the dynamic compression 
factor K, (= E/Eo) is given by 



Q(/C) 
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(15) 



The relations above are valid for gaseous disks. For disks 
that have already been depleted of much of their gas, e.g. ei- 
ther through photo-evaporation by their central star or from 
a nearby O star, the GIs may develop in the dust component 
in a similar way. The stability of a dust particle disk is then 
given by 



7rGEdust 



(16) 



where o"d us t is the velocity dispersion of the dust particles. 
If the disk particles and gas molecules have reached ther- 
modynamical energy equipartition, <7d us t can be estimated 



via the ratio of the dust particle mass m v compared to the 
molecular mass /1 of the gas: 



""dust 



(17) 



For micron-sized or larger particles this is negligible com- 
pared to other causes of velocity dispersion (e.g. gas tur- 
bulences which excite the dust component via gas drag). 
However, it is reasonable to assume the dust component to 
be dynamically much cooler for particles of a certain size 
than the gas component. The dynamical friction damps the 
velocity dispersion of the dust grains over a long time but 
the grains r emain inert to short-time gas perturbations and 
turbulence. iBarriere-Fouchet et alJ <l2005l) show that inter- 
mediate sized dust grains (between about 1 mm and 1 m) 
tend so settle efficiently into a thin dust layer in the mid- 
plane of the disk. If we assume that the disk is already at an 
intermediate age, i.e. a few Myr, so that the dust has partly 
coagulated to at least centimetre-sized grains, a significant 
fraction of the dust will be concentrated in or close to the 
mid-plane with a low velocity dispersion. A partial depletion 
of the gas component will lower the particle size required for 
most efficient settling. Therefore we expect Qdust < Q and 
thus a higher likelihood for GIs within the dust layer once 
the gas has been largely depleted. However, the extension 
of the simulation to this case has a more qualitative than 
quantitative character and needs to be ascertained with hy- 
drodynamical models. The role of the dust is discussed in 
more detail in §|5| 



2.2 Numerical methods 

Because we neglect particle-particle interactions and the 
disk-to-star influence, the gravitational acceleration for each 
test particle depends only on the position of the particle and 
both stars at a given time t, whereas the star-star encounter 
is reduced to a two-body problem. This allows to split the 
motion of the bodies into two parts: 

(i) The hyperbolic (or parabolic) encounter of the stars, 
and 

(ii) the motion of the individual particles in the time- 
dependent gravitational field of both stars. 

For each particle the equations of motion are solved numer- 
ically. In each time step the acceleration is computed from 
the time-dependent stellar positions. 

The time-dependent positions of the two stars are given 
by a (semi-)analytical calculation of the orbits. This is re- 
duced to finding the roots of an ordinary equation. In case of 
parabolic orbits it can be solved analytically, whereas in the 
case of hyperbolic orbits a numerical solution is required. 
For the latter we use a Newton-Raphson method. In order 
to reduce the computational costs, we interpolate the or- 
bits by a cubic spline making use of pre-tabulated orbital 
positions. The equation of motion of the test particles is 
solved by a Bulirsch-Stoer method with adaptive step size. 
To avoid singularities, a gravitational softening with a length 
of 0.001 AU has been applied. 

To obtain the local Toomre parameter the surface den- 
sity is measured by "probe particles" (probes), which sum- 
up within their probe-volume the masses assigned to the 
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test particles. The method is described in more detail in 
Appendix 1X1 



3 SIMULATION OF INDIVIDUAL 
ENCOUNTERS 

To demonstrate the general effects of star-star encounters we 
simulated the fly-by of a star with 0.5 Mq with ro = 120 AU 
pericentre distance and an eccentricity of e = 1.5. For com- 
parison we reran the fly-by with a larger encounter distance 
of r = 150 AU. 

3.1 Initial conditions 

The disk particles are initially set on circular coplanar or- 
bits around the Sun. Depending on the kind of simulation, 
the particle positions are set up as a polar grid arrangement 
(regular distribution) or as a random distribution preserv- 
ing the mean density distribution. For density evaluations 
the regular distribution yields less noisy numerical results. 
On the other hand, the calculation of the mean angular mo- 
mentum (for the estimate of the tilt of the disk) requires the 
radial mass distribution to be as smooth as possible. Thus, 
we used a random reshuffling of the initial particle positions 
for calculating the tilt in § [4.31 emp loying the M ersenne 
Twister algorithm iMatsumoto fc Nishimuralll998f) . 

The initial distance of the Sun and the perturber is set 
larger than 1000 AU to prevent any significant tidal effects 
before the encounter. For the same reason, the simulation 
ends after the distance has reached at least the initial value. 
For higher resolution, the computed population of 250,000 
particles and 1000 probes is limited to the area between 25 
and 30 AU. Such a selection for increasing the spatial reso- 
lution is possible because there are no interactions between 
the particles. The absence of further interactions also allows 
for a simple scaling to other disk radii and densities, so the 
results can be scaled up or down for the orbits of Uranus or 
Pluto and the larger Kuiper belt objects. In a special run, 
which is described in § 15.21 we extend the simulated area to 
100 AU, the known region of the Edgeworth-Kuiper-Belt. 
We also extended the disk area to 50 AU but reduced the 
particle number to 10,000 for the plots shown in Figs. [5] 

3.2 Passage of a 0.5 Mq star at 120 AU 

A close encounter of a 0.5 Mq star at 120 AU with e = 1.5 
gives a good illustration of the fundamental effects of tidal 
perturbations on a circumstellar disk. Fig.|5|shows the over- 
all development of the perturbed disk within 50 AU dur- 
ing and after an encounter at 120 AU periastron distance. 
Within only 400 years strong density fluctuations appear as 
thin tidal arms of highly compressed material. These tidal 
arms persist a few 10 2 years up to 10 3 years, depending 
on the mass and the encounter distance of the perturber. 
Within the spiral arms, material can be carried far outwards 
and even become unbound. The time-dependent surface den- 
sity is shown in Fig. [3] for encounter distances of 120 AU 
and, for comparison, 150 AU. It can be seen that the fly-by 
distance has a strong influence on the magnitude of the E 
peaks; while a 120 AU fly-by causes an increase of more than 
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Figure 2. Disk perturbations during hyperbolic fly-by (e = 1.5) 
of a star with 0.5 Mq at ro = 120 AU. Only the regions within 
an initial radius of 50 AU are shown. 

a factor of 10, the 150 AU encounter does little more than 
doubling E. 

As the fly-by proceeds we follow Q in each counting 
volume and at each time step the minimum Q is selected. 
A plot of the Qmin values within the computed disk area is 
shown in Fig. 2] for the isothermal (adiabatic index 7=1) 
and the IHh/He adiabatic case (7 = 1.4). A disk with a lo- 
cally isothermal equation of state and an initial Q — 8.5 ex- 
hibits regions with Q < 0.8 within less than 500 yr after a 
perihelion passage closer than about 120 AU, satisfying the 
Toomre instability criterion (Q < 1). For the adiabatic bor- 
der case, Q locally falls to 1.2, still below t he upper lim it for 
Toomre instability Q « 1.7. According to iBossl (2004), the 
combined convective and radiative cooling is sufficient to 
form clumps similarly to the isothermal case. Cooling times 
are of the order the orbital period even for optically thick 
disks an d there ful f il the criterion for fragmentation investi- 
gated by iGammid (2001) . Thus, an isothermal equation of 
state is expected to be a reasonable approximation. 



4 RESULTS OF THE PARAMETER STUDY 

We surveyed the parameter space (m, e, ro) in order to find 
the regions where Q < 1 locally occurs. The parameters have 
been varied successively with appropriate step sizes (Tab.0. 

Here we describe the major results of our parameter 
study dealing with the triggering of gravitational instabili- 
ties, the masses and orbits of seed candidates for coplanar 
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Figure 3. Surface density E vs time for ro = 120 and 150 AU 
measured in one probe at a radius of 30 AU. In the closer fly- 
by E reaches easily 10 times the initial value about 500 yr after 
periastron. The wider fly-by leads only to slight density variations, 
which illustrates the strong dependence on ro- The peaks in E are 
due to the density waves moving with a different velocity than the 
Keplerian orbital motion. 
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Figure 4. Minimum values of Q within the region 
25 AU < r < 30 AU (initial conditions) during a 120 AU 
encounter of a 0.5 Mq perturber for the isothermal and the 
adiabatic case (7 = 1.0 and 7 = 1.4, respectively). Even in 
the adiabatic case Q grazes the zone of Toomre instability 
(Qmin ~ 1-7), whereas the more realistic isothermal case almost 
ensures the Toomre instability criterion (Q m i n < 1). 

prograde passages, and the angular momentum tilt due to 
inclined passages. The possibility of triggering for inclined 
passages has been studied for isolated cases. 

4.1 Triggered Toomre instability for prograde 
passages 

Our parameter study for prograde coplanar passages shows 
that the Toomre instability criterion Q < 1 is locally satis- 
fied for all stellar masses above 0.1 Mq and all eccentricities 
between 1 and 10 for encounter distances ro < 80 AU. For 
perturber stars with 1 Mq and e = 10 we found Q < 1 
occurrence even for fly-bys with ro < 170 AU. This critical 
encounter radius ro, C rit, the upper limit for ro up to which the 
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Table 1. Ranges and step sizes for the varied parameters m, e 
and ro 
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3 5 


10 


50 


0.1 


0.77 


0.71 


0.61 


0.53 


0.44 


0.26 


0.3 


0.73 


0.68 


0.58 


0.50 


0.41 


0.25 


0.5 


0.69 


0.64 


0.55 


0.48 


0.39 


0.24 


0.7 


0.67 


0.62 


0.53 


0.46 


0.38 


0.23 


1.0 


0.63 


0.58 


0.50 


0.44 


0.36 


0.21 



Table 2. Coupling radius ag in units of ro, depending on the 
mass and orbit of the perturber 

Toomre instability criterion is satisfied, depends on both the 
perturber mass and the encounter eccentricity in a mono- 
tonically increasing manner, as can be seen in Fig. Q3 This 
reveals two major trends: Increasing the mass leads to an 
increased maximum perihelion distance ro,crit, and also a 
higher eccentricity of the perturber orbit results in larger 

^O.crit ■ 

The mass dependency is simply caused by an enhanced 
perturbation with increasing mass. The eccentricity depen- 
dence on the other hand is less obvious. Unless the encounter 
velocity is extremely high or the periastron distance ex- 
tremely small (i.e. close to the particle orbits) the angular 
velocity of the perturber is normally significantly smaller 
than that of the disk particles. Thus an increment of the en- 
counter velocity for a given periastron radius yields a larger e 
and a smaller particle-perturber relative velocity, i.e. longer 
exposure time. 

An important question is how much mass is included 




Figure 5. Upper border ro cr it of the regime of Q < 1 occurrence 
in the (m, e,ro) parameter space for the isothermal equation of 
state. Q < 1 occurs for values of ro below the border indicated 
by the grid surface. Clearly visible is the increasing ro trend for 
increasing perturber mass and eccentricity. The latter trend is 
caused by the decreasing relative angular velocity between per- 
turber and test particles. The key on the right indicates the ro- 
contour levels plotted in the m, e-plane. 
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Figure 6. Critical upper pericentre distance ro jC rit for a sam- 
ple of encounters with non-zero inclination i compared to the 
corresponding coplanar prograde fly-by. The perturber mass is 
0.5 MQin all cases. Retrograde passages result in ro jCr it about 
50 % of that for zero inclination independently of the eccentricity. 
jp is the angle between the pericentre and ascending node. 



in a typical probe with observed Q < 1. A typical probe in 
our simulation has a radius h = 0.1 AU (Appendix [XJ, but 
the total region of a local Toomre instability could be larger 
than that. In an unperturbed disk E has a value of about 
40 gem -2 or 1.5 M® AU -2 at 30 AU, which results in about 
0.03 Earth masses (M$) within h. Within a tidal arm E can 
reach 10 to 20 times this value, yielding 0.3 to 0.6 per 
probe. 

Since the probe contains only a fraction of the Toomre- 
unstable material of the tidal arm, the total mass of the 
unstable region is likely to be much larger. A typical com- 
pressed arm is about 0.5 AU wide. A spherical region of 
0.5 AU contains 10-20 M® or about one Neptune mass. If 
one per cent of this is dust dNattall200j) a body of about 
0.1-0.2 M® can form out of this seed. This would skip much 
of the coagulation time otherwise needed, and is about 100 
times more massive than the most massive KBOs. 



4.2 Triggered Toomre instability for retrograde 
and inclined passages 

In a few examples we also studied the triggering effects of 
an inclined prograde and a retrograde passage of a 0.5 Mq 
perturber. We varied the inclinations between 45° and 180° 
in steps of 45°. As one would expect, the perturbation in 
both cases is smaller than that of a coplanar prograde pas- 
sage being due to the larger distance between the perturber 
and that edge of the disk with the lowest angular velocity 
relative to the perturber. The results are summarised in Fig. 

For the coplanar retrograde passage (inclination 
= 180°), we found ro, C rit = 65 AU for e = 1.5 and 
^o.crit = 75 AU for e = 10. Both are about half the ro, cr it in 
the coplanar prograde case, i.e. the reduction of ro.crit due to 
inclination is apparently independent of the eccentricity. For 
other inclinations the decrease of ro, C rit depends on the an- 
gle tp between the pericentre and the ascending node, where 
ip — 90° leads to significantly higher ro.crit than ip = 0. 



t- -2500 vr 
d = 1 094 A1J 



t - yr 
d = 1 00 AU 




■f= 500 yr 
d - 334 AU 



t- -100 vr 
d= 125 AU 




H 1 h 



t= 5000 yr 
d - 1 787 AU' 



-200 -100 100 200 -100 100 200 

x I AU x / AU 

Figure 7. Parabolic fly-by of a 0.5 Mq star at 100 AU with 45° 
inclination. The angle tp between pericentre and ascending node 
is 90° . 



4.3 Tilt of the angular momentum vector 

In inclined encounters there is the possibility for tilting the 
plane of the disk due to the transfer of angular momentum to 
the orbit of the perturber. We investigate the possible origin 
of the obliquity of the Sun's rotational axis relative to the 
ecliptic for inclinations i between 5° and 175°and ip = 0. 
Passages with inclination less than 90° are designated as 
prograde, and those beyond 90° as retrograde. In this section 
the radius of the model disk is set to 100 AU to measure the 
angular momentum change of all parts of the disk, because 
we expect the resulting tilt to be close to the mean tilt of 
the particle orbits and a global rather than a local effect. 
The global effects on the disk are shown in Fig. Q 

After the encounter, when the distance between the 
stars is large enough to prevent significant further tidal ef- 
fects, the sum of the orbital angular momenta of all disk 
particles is calculated. The results are plotted in Fig. |H| the 
most remarkable feature is that the tilt has its maximum of 
13° for retrograde passages of about 140° and ro = 75 AU, 
whereas the local maximum in the prograde area at 45° ap- 
pears much lower for close encounters below 100 AU (6.5° 
for 75 AU). For less close encounters, however, the maximum 
tilts are essentially equal for both prograde and retrograde 
passages. In both cases, an encounter below about 125 AU 
results in a maximum tilt of more than 6° which is close to 
the observed solar tilt of 7°. We return to this point in § 15.31 
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15 30 45 60 75 90 105 120 135 150 165 180 

i/° 

Figure 8. Tilt /3 of the disk after fly-by of a 0.5 Mq perturber as 
a function of inclination and closest distance, i < 90° are prograde 
encounters, i > 90° retrograde ones. 



5 DISCUSSION 

5.1 Cooling and fragmentation of disks 

Based on an analytical estimate lRafikov| i2005t 



found that 

the fragmentation of protoplanetary disks requires unusual 
high surface densities and temperatures. According to his 
paper a local Toomre instability does not necessarily lead 
to fragmentation. It is crucial that the instability can be 
achieved despite temperatures which are high enough to pro- 
vide efficient radiative coolin g via the Stefan-Boltzmann law. 
According to Rafikg3 <l2005f) the temperature T cr it required 
for efficient cooling at 10 AU ranges between 200 and more 
than 1000 K depending on the disk model. On ly if the disk 
mass is a factor of 10 higher than that used bv lBoss! (I2004T) 
is a Toomre instability possible. Because, even if the initial 
temperature of the condensed region is less than T cr i t , the 
cloud being more unstable, the rapid increase of temperature 
during gravitational contraction eventually compensates the 
increasing pressure. 

The minimum disk surface density E m j n a nd the critical 
temperature T cr it are given by 



E > E min = fl 7/5 {nGQo 



and 



-6/5 



1 

0? 



-, 1/5 



/(r)] 1/5 (18) 



T > Tent = ft 4/5 [C-kQoGg) 



-2/5 



3/5 



lf(r)} 



2/5 



(19) 



where Qo w 1 is the upper Toomre limit for disk instability, 
£ = 2£(7 — 1) is a parameter to be satisfied by the cooling 
time f C ooi and the orbital fr equency O, Qt c oo i < £ ~ 1. Note 
that our T cr it equals T m i n in Rafikovl i2005f) . a is the Stefan- 
Boltzmann constant, k the Boltzmann constant, /i « 2.3 u 
the molecular mass of the disk gas and 



fir) 



T + 



(20) 



where r is the optical depth which we estimate to be in the 
order of several tens in a compressed region. The value of r 
depends on the opacity and the surface density of the par- 
ticular disk region. Given the model parameters used in this 



paper, an estimate of r = 60 of the compressed disk region 
and £ = 1, conditions 1181 and I19H require E > 630 gem -2 
with T cr it = 330 K at 30 A U. Numerical simulations sug- 
gest £ w 3 iGammiel 1200 J) ; in that case the limits are 
Emin = 510 g cm~ 2 and T cr it = 210 K. Assuming a tidal com- 
pression by a factor of 15 (which is a reasonable value in our 
simulations) the surface density condition is easily achieved 
in the compressed area, even with sig nificantly lower dis k 
masses than proposed bv lBoss! (I2004T) or lMaver et alJ i2004) . 
The assumed optical thickness in these cases corresponds to 
a mean opacity lu = r/E « 0.1 cm 2 g" 1 of the disk mate- 
rial, which is a reaso nable value for condensed region s in 
protoplanetary disks JPodolakll20ol iPollack et alJll985l for 
centimetre-sized dust grains). 

Another problem which we do not treat in this paper 
is the transformation of the presumably unstable regions 
into (proto-)planets. In case of the formation of Neptune or 
Uranus this includes the question of how to get rid of a large 
amount of gas with respect to the dust, by this ending up 
with the observed high mass fraction in solid material in 
these planets. 

In principle, the following scenarios for the evolution of 
a gas-rich dust disk are possible: 

A dust phase can partly decouple from the gas, re- 
sulting in a comparatively thin dust layer jNattal 120041 : 
iBarriere-Fouchet et al J 12005). Since the dust evolves basi- 
cally pressure-free it is much more prone to gravitational 
instability than the gas and, hence, instabilities in the dust 
phase could grow faster and form denser regions than the 
corresponding instabilities in the dynamically hotter gas. A 
plausible way of increasing the global solid-to-gas ratio is 
as follows: Jupiter and Saturn form traditionally within a 
few Myr. Meanwhile the disk is partially photo-evaporated 
such that the dust-to-gas ratio increases significantly, at 
which time the fly-by-perturbation occurs. From Fig. |S] 
we see that gravitational instabilities can be induced for 
fo.crit < 200 AU, and from Fig. it can be seen that the 
probability of such an encounter lies between 10 and 60 % if 
the Sun was born in a cluster with a radius of ~0.5 pc con- 
taining 1000 < N < 10000 low-mass stars (one encounter 
within ~10 Myr). In this scenario the GIs would be induced 
in a dust-dominated disk (see § 12. H . 

A preferential loss of gas locally might also be possi- 
ble when the clumps survive sufficiently long and mass seg- 
regation by gravity brings the dust to the center of the 
clumps, by this shielding the dust from a later loss by 
photo-evaporation. Boss (2004) suggested that the clumps 
can survive for a few dynamical periods which might 
be sufficient for the required mass segregation. Moreover, 
mass segregation can be pos sibly sped-up by vortices, as 
iKlahr fc BodenheimeJ J2jp04) show. The solid component 
tends to settle r apidly in the center of a vort ex due to dy- 
namical friction. IKlahr fc BodenheimeJ <|2003T) describe the 
formation of vortices from baroclinic instabilities. However, 
we expect vortices to occur also in gravitational instabilities 
if there is an initial spin angular momentum. The contract- 
ing GI will then be spun-up by the Coriolis forces. 

In a border case the gas fraction may have been lowered 
due to normal depletion mechanisms during disk evolution. 
In the case that the dust may have already settled in or near 
the mid-plane, the horizontal centring due to vortices will 
result in an extreme concentration of dust within a small 
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volume, thus skipping the otherwise necessary vertical con- 
traction due to self-gravity. 

We want to stress that the discussed possibilities for 
getting high dust mass fractions in the clumps (i.e. an ini- 
tial value problem or an evolutionary process) cannot be 
decided on the basis of our simple simulations. A more re- 
fined multi-component treatment of dust and g as as separate 
phases (e.g. as done bv lTheis fc Orloval ll2004) in the case of 
dusty galactic disks) would allow for investigating these in 
detail, but this would be part of a future investigation. We 
do believe, however, that our found unstable regions would 
be excellent large-scale seeds acting as birth places of (proto- 
) planets. 



5.2 Orbits of candidate planetesimals 

On a more speculative note our simulations allow us to study 
the likely initial shapes of the orbits of candidate protoplane- 
tary seeds, which are chosen here as the probes where Q < 1, 
indicating the possible formation of a planetary core. The 
flyby distance is 100 AU, and two cases, a coplanar and an 
inclined orbit (i = 45°) are investigated. The perturbations 
have their greatest impact on the region beyond 30 AU, i.e. 
the Kuiper belt region. 

The distribution of the eccentricities e and semimajor 
axes a is shown in Fig. [0] where from 10000 test volumes ini- 
tially at radii between 20 and 100 AU only those are plotted, 
for which Q < 1 occurred at least once. 

It can be seen that the orbit of the KBO Sedna is not 
far outside the bulk of candidate seeds, and well contained 
in a less-dense peripheral region of candidate seeds for the 
coplanar case. It is still close to it in an 45° encounter. Also 
Pluto is contained in this population, despite its resonance 
with Neptune, which most probably is the result of orbital 
evolution after the formation. Neptune, however, lies below 
the candidate seed population but close to a trail of low- 
eccentric particles. Also due to its low eccentricity Quaoar 
lies outside the population. The eccentricities in this plot 
range from 0.02 (nearly circular) to 1 (parabolic ejection) 
and would even exceed 1 if plotted over the perihelion dis- 
tance rather than the semimajor axis. Even the recently an- 
nounced large KBO 2003 UB313 fits this plot. However, its 
the longitude of its ascending node deviates more from those 
of the other shown objects. Note that, if the present propo- 
sition is realistic, we expect our solar system to be void of 
KBOs below and a bove the region p opulated in Fig. 

We note that llda et alJ i2000f) propose that the eccen- 
tric orbits of the Kuiper-Belt objects may have formed as 
a result of tidal deformation of a pre-existing planetesimal 
disk. Our alternative notion is that the eccentric outer-solar- 
system objects may be a result of an encounter which lead 
to the formation of Uranus, Neptune and the KBOs. 



5.3 Earlier SPH simulations of tilted disks 

A more extensive study of the disk tilt due to stellar encoun- 
ters than ours (§ I4.3P has been made by Heller (1993), who 
used an SPH model with disk-star and disk-perturber inter- 
action instead of restricted iV-body computation. His results 
are essentially the same as in our study. Minor differences 
due to the absence of viscosity are possible. According to 
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Figure 9. Distribution of semimajor axes a and eccentricities e 
of the Q < 1 probes (indicating candidate planetary seeds) after 
a parabolic fly-by at 100 AU with an inclination of i = 0° (upper 
panel) and i = 45° (lower panel). The orbits of some actual major 
bodies are included as well. Pluto, Sedna and 2003 UB313 are 
contained in or close to the main population, whereas Quaoar 
lies outside and Neptune below the "trail" at 25-30 AU. If the 
early Neptune had a slightly more eccentric orbit (e = 0.07) it 
would have fit this plot. 



Heller (1993), the ability of the disk to deform, which is 
maximal for non-interacting test particles, enhances the re- 
sulting disk tilt. The perturbation of the perturber-star mo- 
tion due to disk-star interactions also alt ers the result s, but 
since the effects are small in most cases (lHelleJll993f) . they 
can be neglected here. The important factor causing Heller's 
results on the tilt to be somewhat larger than ours is due to 
his large particle removal radius of 800 AU, whereas we ne- 
glected particles with a semimajor axis larger than 150 AU. 



5.4 Limits and perspectives 

The restricted iV-body method has its great advantage in 
low computational efforts. This allows us to perform a wide 
parameter study within an acceptable computing time, but 
on the other hand this also limits the accuracy of the an- 
swers. 

Since we do not model self-gravity our approach is only 
reasonable if the mass of the disk is small compared to the 
solar mass, making the latter the dominant source of gravity. 
The same applies for the dynamical feedback on the orbit 
of the stars and the wobbling of the Sun (which would alter 
the density evolution of the disk) due to asymmetric disk 
perturbations. 

The absence of hydrodynamical interactions limits the 
simulation time within which reasonable results for the den- 
sity evolution can be expected. Furthermore, the effect of 
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coolin g mechanisms as mentioned bv lBossI i2002T) and iBossl 

(2004) can only be approximated by adjusting the adiabatic 
exponent 7. 

Due to these restrictions this work provides only a rough 
estimate of the expected magnitude of triggered instabili- 
ties. The evolution of the strong initial perturbations is de- 
termined mainly kinematically. So the results are reliable if 
only the immediate effects are of interest. As time proceeds 
the evolution will be more and more determined by hydro- 
dynamical and other interactions within the disk. But this 
uncertainty mainly affects the borders of the candidate pa- 
rameter subset rather than the gravitational-triggering sce- 
nario as a whole. 

Despite these limitations, our study gives us an impor- 
tant guideline for a more tightly focussed but computation- 
ally extremely expensive hydrodynamic research. 



6 CONCLUSIONS 

Using a simple, but fast numerical method a parameter 
study has been made to test the possibility of gravitation- 
ally triggered planet formation. This scenario is supported 
by analytical estimates of the relatively high encounter prob- 
ability in a typical young cluster (Fig. ^ as well as by the 
observed obliquity of the solar rotational axis being a possi- 
ble consequence of an inclined fly-by. 

As the main outcome of this study we summarise the 
following results: 

(i) Prograde passages of stars with masses of 0.1-1 Mq 
can trigger Toomre instabilities within a 30 AU radius in a 
disk of 0.07 Mq if an approximately isothermal behaviour 
of compression is assumed. Toomre instabilities are also ex- 
pected for adiabatic compression. Higher perturber masses 
and eccentricities increase the density perturbations and 
therefore allow for larger areas to become unstable. 

(ii) Retrograde passages cause local Q < 1 for smaller ro 
(about half the value of the prograde case). Inclined passages 
yield Q < 1 for intermediately close fiy-bys. 

(iii) A fly-by of a 0.5 Mq star with an inclination between 
30° and 60° and ro < 100 AU tilts the disk plane by 3°-6°. 
This is very similar to the observed solar obliquity. 

The duration of the locally triggered gravitational in- 
stability is in the range of centuries to millennia. There- 
fore, if these Toomre instabilities can lead to fragmentation 
and enhanced core-formation they may solve the discrep- 
ancy between the short disk lifetime and the long (classical) 
formation times of the outer planets and maybe even the 
Kuiper-Belt objects (KBOs). Even if only a fraction of to- 
day's Neptune mass was assembled by triggering, this might 
be the required seed for the ongoing coagulation and ac- 
cretion to be finished before the disk dissolves. The high 
solid-to-gas ratio of Uranus and Neptune may be explained 
in this scenario if their formation was induced by a fly-by 
at a solar-system age > 5 Myr when the disk was depleted 
of gas. In this situation Qdust < Q such that the dust disk 
should be even more unstable towards GIs than suggested 
above. This sets possible constraints on the physical prop- 
erties of the Sun's birth cluster. 

Our assumed disk mass is comparable to the minimum 
solar nebula mass. This avoids implausible assumptions on 



the disk mass as they are required for pure gravitational 
fragmentation of unperturbed disks. 
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APPENDIX A: DENSITY EVALUATION 

Calculation of the local Toomre disk stability requires the 
evaluation of the local surface density. At first, all disk par- 
ticles are weighted by an individual mass mj which corre- 
sponds to its initial orbit radius rt, in such a way that 

JV 

™b = M d (n < r < r 2 ) , (Al) 

i 

where r\ is the inner and r 2 the outer border of the simulated 
disk region. Since the N particles are initially positioned on 
a polar grid consisting of n r rings with n v particles each, 
which results in a r _1 law of the surface density E, the mass- 
weighting must be proportional to r -1 / 2 to fit the needed 
r" 3/2 law (EqE2J. Together with eqn. llA"Tt this leads to 

m(r b ) = ^^- ^"^L -Md, (A2) 

-1/2 



E' V -1/2 V^n 



r 

6=1 



where r m i n and r max are the inner and outer radius of the 
real disk. Here we use r m i n = 3AU and r max = 100 AU. 

The density evaluation is achieved by simply adding the 
nit within "counting-volumes" surrounding massless probe 
particles which move along with the disk particles on ballis- 
tic (Keplerian) orbits. Within the counting-volume the den- 
sity is calculated in SPH-style using 

JV 

{Q) w =^m h W(d b ,h), ( A3 ) 

6=1 

where d b is the distance of the b-th mass-weighted test par- 
ticle from the center of the probe particle. The kernel func- 
tion W provides a smooth distance-weighting of the test 
particles for reduction of numerical noise. In our model t he 
compact spline kernel from lMonaghan fc Lattanziol ill 9851) is 
used. The choice of the kernel radius h determines the mean 
number of particles nh within the probe volume. For a ran- 
dom distribution, nu > 100 is required to keep the noise low 
enough. A regular distribution, as we used here, shows ap- 
propriate results even for ~ 10. For a given (average) 
and disk population parameters ri, r 2 and N we get h from 



n h (ri + r 2 ) (r 2 - n) 

N ■ (A4) 

Thus, typically, h = 0.1 AU for n = 25 AU, r 2 = 30 AU and 
N = 250000. 

Since we want to identify unstable regions, i.e. regions 
where Q < 1 locally, we applied about 1000 of these test 
volumes, which concentrate at the areas of the highest den- 
sity. Their dynamic behaviour is the same as that of the disk 
test particles. Thus, a high spatial resolution at the areas of 
interest is automatically provided. 



